34 ,h AIAA Fluid Dynamics Conference 


AIAA 2004-2539 


Numerical Simulation of Roughness-Induced 
Transient Growth in a Laminar Boundary Layer 




* 


Paul Fischer* 

Argonne National Laboratory \ Argonne , JL 60439 
and 

Meelan Choudhari f 

NASA Langley Research Center Hampton , VA 2368 1 


Numerical simulations are used to examine the roughness-induced transient growth in a 
laminar boundary-layer flow. Based on the spectral element method, these simulations 
model the stationary disturbance field associated with a nonsmooth roughness geometry, 
such as the spanwise periodic array of circular disks used by White and co-workers during a 
series of wind tunnel experiments at Case Western Reserve University. Besides capturing 
the major trends from the recent measurements by White and Ergin, the simulations 
provide additional information concerning the relative accuracy of the experimental findings 
derived from two separate wall-finding procedures. The paper also explores the dependence 
of transient growth on geometric characteristics of the roughness distribution, including the 
height and planform shape of the roughness element and the ratio of roughness size to 
spacing between an adjacent pair of elements. Results are used for a preliminary assessment 
of the differences between recently reported theoretical results of Tumin and Reshotko and 
the measurements by White and Ergin. 
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roughness height 

Measure of streamwise velocity perturbation associated with any Fourier mode involved in 
transient growth 

constants involved in scaling of transiently growing perturbation amplitudes with roughness height 
integrated energy for spanwise mode of harmonic index k at a given streamwise location (= \U k 1 2 drj) 
height of roughness element 

ratio of spanwise wave number corresponding to Fourier mode of interest to fundamental spanwise 
wave number based on element spacing within roughness array 
order of discretization polynomial within each element 
power spectral density 

Reynolds number based on roughness height and mean-flow velocity at this height 
Reynolds number based on U «*, and x 

radial distance from center of roughness element, i.e. ((x- Xo) 2 +z 2 ) I/2 
roughness radius 
root mean square 

= Modal profiles of streamwise, wall-normal, and spanwise velocity perturbations corresponding to 
spanwise mode of wavelength X - Xo/k 
free-stream speed 

velocity components along streamwise, wall-normal, and spanwise directions 

streamwise, wall-normal, and spanwise coordinates (x = 0 corresponds to plate leading edge, i.e., 

the virtual origin of the boundary layer in an experiment; z = 0 is aligned with the center of the 
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*o 

Xexp 

Po 

A 


roughness element) 

roughness location relative to plate leading edge 

streamwise coordinate relative to physical leading edge in the experiment of White and Ergin (2003) 
(x exp = x + 70 mm) 

Wave number of the Ao = Aomode, non-dimensionalized by the local similarity length scale x Re x ' 1/2 
Blasius similarity variable (y/x Re x 1/2 ) 
span wise wavelength of Fourier mode 

spacing between adjacent pair of elements from spanwise roughness array 


I. Introduction 

T HE transient growth phenomenon refers to an algebraic amplification of small-amplitude disturbances prior to 
an exponential decay farther downstream. Physically, transient growth arises due to the “lift-up” effect 
associated with streamwise vortex motions (Landahl 1980). Mathematically, transient growth is a direct 
consequence of non-normality of the linear stability operator (Trefethen et al. 1993). Nonnormality implies that a 
suitable superposition of linearly stable eigenmodes may involve significant cancellation between the associated 
flow perturbations. Due to streamwise variations in the underlying mean flow (i.e., nonparallel effects), the 
destructive interference characteristics cannot be sustained while the above set of disturbances evolves in the 
downstream direction. The gradual weakening of mutual cancellation results in a transient growth in perturbation 
amplitudes before the exponential decay of the stable eigenmodes eventually comes into play. 

For suitable inflow conditions, the initial cancellation is strong enough for the transient amplification ratios to 
meet or even exceed the range of growth factors that is known to correlate with the onset of transition due to the 
exponential growth of linearly unstable eigenmodes. Consequently, the transient growth paradigm has emerged as 
an alternate scenario for laminar-turbulent transition, especially under subcritical (i.e., linearly stable) flow 
conditions (Reshotko, 2001). An upper bound on the transient growth ratios is provided by the optimal growth 
theory (Farrell 1988, Butler and Farrell 1992, Andersson et al. 1999, Luchini 2000, Tumin and Reshotko 2001). 
According to this theory, the optimal initial conditions are associated with purely wall normal and spanwise velocity 
perturbations resembling a spanwise array of streamwise vortices; the transient growth occurs as these initial 
conditions evolve into streak like motions that are dominated by perturbations in the streamwise velocity. 

The physical relevance of optimal growth theory to boundary-layer transition is largely dependent on receptivity 
characteristics of the laminar boundary layer, in particular, whether and how such optimal (or near-optimal) initial 
conditions can be realized in a natural disturbance environment. Transient growth of low-frequency boundary-layer 
disturbances due to free-stream turbulence of weak through moderate intensity has been documented in a number of 
experiments (Kendall 1985, Westin et al. 1994) and also predicted using theoretical models (Andersson et al. 1999, 
Leib et al. 1999). The work by Leib et al. (along with the related earlier works by Choudhari (1996) and Bertolotti 
(1997)) illustrates the need to model the transient growth phenomenon as an inhomogeneous boundary value 
problem associated with a physically realizable forcing environment, rather than via the optimal growth formulation 
alone. 

The analogous issue of whether initial conditions appropriate for near-optimal growth can also be seeded via 3D 
surface roughness (and, if so, of what kind) has received less attention thus far. Reshotko and Tumin (2002) 
recently used a linear amplitude criterion based on optimal disturbance growth to correlate with an existing database 
for subcritical transition on rough axisymmetric nosetips. The finding that this linear criterion was able to correlate 
the major trends involving the effects of roughness height and surface temperature reinforced the notion that surface 
roughness might also be an important catalyst for inducing transient growth phenomena. However, because the 
maximum roughness heights associated with the experimental database were comparable to the boundary-layer 
thickness, the reasons behind the success of a transition criterion based on a purely linear model are not clear at this 
stage. 

Quantitative measurements of the stationary disturbance field induced by 3D roughness elements in the form of a 
circular disk were carried out by Tani et al. (1962), and then by Kendall (1982) in conjunction with his pioneering 
measurements of the transient growth due to weak free-stream turbulence. Due to the difficulty in measuring the 
weak velocity perturbations induced by an isolated roughness element, Kendall had to use signal-averaging 
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techniques to enhance the accuracy of the measured data. Gaster et al. (1994) circumvented these challenges by 
using a surface diaphragm oscillating al very low frequency. Like Kendall, Gaster et al. observed a switchover in 
the spanwise disturbance profile from a unimodal behavior immediately behind the roughness element to a profile 
involving multiple spanwise extrema at farther downstream locations. However, the respective datasets are found to 
differ with respect to the number and the signs of these spanwise peaks. 

Gaster et al. noted that the appearance of additional spanwise peaks was accompanied by a significant transient 
streamwise growth within a band of low wave number disturbances (compared with the peak wave number of the 
unimodal disturbance profile). The transient growth phenomenon was also captured in the direct numerical 
simulations by Joslin and Grosch (1995), who modeled the low-frequency oscillation of the diaphragm from Gaster 
et al. as a stationary hump. However, an uncertainty in the mode shape of the vibrating diaphragm in the experiment 
prevented an in-depth quantitative comparison with the measured data. 

The investigations by Kendall (1982), Gaster et al. (1994), and Joslin and Grosch (1995) were carried out prior 
to the development of a spatial theory for optimal growth (Andersson et al. 1999, Luchini 2000, Tumin and 
Reshotko 2001). More recent wind tunnel experiments by White (2002) and White and Reshotko (2002) also 
demonstrated the occurrence of roughness-induced transient growth in the intermediate wake region for a variety of 
roughness configurations involving cylindrical (i.e., circular disk) roughness elements and sand-grain roughness. 
Although the measurements of White (2002) were consistent with some aspects of the optimal growth theory, the 
magnitude of the transient growth associated with a spanwise periodic array of cylindrical elements was found to be 
suboptimal. In reducing the hot wire data to determine the spanwise variation in the mean velocity. White (2002) 
used a linear extrapolation technique to estimate the local plate surface location at a given spanwise coordinate. 
Subsequent analysis by White and Ergin (2003) showed that the distorted mean-flow profile near the wall has a 
finite curvature within the region of transient growth and, therefore, some of the results reported in White (2002) 
may not have been correct. White and Ergin presented additional data for a spanwise array of cylindrical roughness 
elements that was reduced by using a spanwise global method for estimating the plate surface location. These new 
measurements further quantified the variation in the amplitude of algebraically growing disturbances with respect to 
the roughness height. 

Tumin and Reshotko (2004) have presented a clever theoretical analysis of roughness-induced transient growth; 
they used a generalization of the method that was developed and applied earlier for receptivity predictions in the 
context of linearly unstable Tollmien-Schlichting modes. Using biorthogonal eigenfunctions of the linear, quasi- 
parallel disturbance equations, they were able to predict the stationary disturbances in the wake of an array of 
roughness elements without having to solve for the disturbance field either above or upstream of the roughness 
array. The disturbance profiles obtained with this procedure were used as initial conditions for a marching 
procedure to confirm the realizability of optimal growth for this configuration. Even though the flow conditions 
used in their case study were nominally identical to those of the White and Ergin (2003) experiment, Tumin and 
Reshotko’s results show some differences in the characteristics of the observed transient growth. It is not clear 
whether or to what degree those differences may be attributed to the differences in planform shape (namely, square 
roughness elements in the analysis versus circular ones in the experiment) or the departure in experiment from 
theoretical assumptions (linearity, and parallel and quasi-parallel flow assumptions in the nearfield and farfield, 
respectively). 

In summary, the qualitative features of roughness-induced transient growth are somewhat consistent between 
theory and experiments. However, additional theoretical studies are necessary to establish quantitative agreement 
for a common set of flow and roughness parameters. Such studies are also desirable in view of the various 
challenges in making reliable measurements of small-to-moderate amplitude spanwise variations in the mean 
velocity field, so as to help resolve some of the discrepancies between previously reported experiments involving 
roughness-induced disturbances. To help address these issues, this paper presents numerical simulations for a select 
set of test conditions from the experiments of White and Ergin (2003). A secondary objective behind this effort is to 
extend the significant body of previous numerical investigations from roughness elements with a smooth shape (i.e., 
C° and smoother) to the cylindrical, disk-shaped roughness elements that are encountered in a number of research 
studies as well as technological applications. 

An outline of the paper is as follows. In section II, we outline the flow configuration and the numerical 
algorithm used for the simulations. Section III describes the results obtained for the baseline configuration, with an 
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emphasis on detailed comparison with the White and Ergin (2003) data as well as on results pertaining to How 
quantities that cannot be easily measured in an experiment. Findings from a rather limited parameter study 
involving additional roughness configurations are outlined in section IV. Concluding remarks are presented in 
section V. 


II. Flow Configuration and Computational Procedure 

The primary flow configuration used for the numerical simulations corresponds to a flat-plate boundary layer 
with a spanwise periodic array of circular-disk roughness elements located at the distance of x 0 = 230 mm from the 
sharp leading edge. In comparison with free-stream turbulence, which provides sustained external forcing 
throughout the streamwise extent of the flow, a localized roughness array is perhaps more relevant to the initial 
value problem modeled by the optimal growth theory. The diameter of the roughness elements is equal to 6.35 mm, 
with the spacing between each pair of adjacent elements being Ao = 19 mm, i.e., approximately 3 times the 
roughness diameter. The free-stream speed is lOm/s. The baseline simulations correspond to a roughness height of 
h = 0.57mm (i.e., Re k = 119). These values correspond to case 3 in the experiments of White and Ergin (2003), 
which was the focus of most of the results presented in their paper. 

The Nek5000 code developed at the Argonne National Laboratory (Fischer et al. 2002) was used for the 
numerical simulations. In this code, the unsteady incompressible Navier-Stokes equations are integrated in time by 
using the spectral element discretization in space and a third-order, operator-splitting formulation in time. The time 
advancement scheme includes an implicit, linear Stokes solver combined with an explicit subintegration of the 
nonlinear advection term. The discretized Stokes problem is solved by using operator splitting, which results in 
decoupled equations for the velocity components and pressure. Despite the fact that only C° continuity is imposed 
across element interfaces, the overall results are known to be spectrally accurate in space. The combination of 
geometric flexibility and spectral accuracy makes the spectral element method particularly well-suited for simulating 
the present flow configuration, as well as to meet the long-term goal of transition simulations involving complex 
geometry such as multielement, high-lift configurations. 

To facilitate a series of simulations for multiple roughness heights, the computational domain is split into a 
nearfield region including the roughness array and a farfield region occupying the smooth downstream section. The 
nearfield subdomain is discretized via an unstructured array of hexahedral elements, some of which have curved 
boundary faces, while the farfield domain is composed of a regular array of hexahedral elements that admits the use 
of fast tensor-product solvers. The interface between the nearfield and farfield zones is a conforming one, so that the 
adjacent elements on opposite sides of the interface align perfectly with each other. Figure 1 shows a schematic of 
the grid topology, including a close-up of the grid in the vicinity of the roughness array. The nearfield and farfield 
portions of the grid contain a minimum of 762 and 2160 elements, respectively. The nearfield domain in the 
baseline case extends from 19 r 0 upstream of the roughness location to 15 r 0 in the downstream direction, whereas its 
free-stream boundary is located at 20 r 0 , i.e., nearly 21 times the boundary-layer thickness at the roughness location. 
The farfield domain begins at 8 mm ( ~ 2.52 r 0 ) downstream of the roughness location (x 0 ), and extends up to x « x 0 
+ 500 mm. Progressive stretching of the wall-normal grid with increasing x provides approximately constant 
resolution across the boundary layer at all streamwise locations. The streamwise grid is stretched in the downstream 
direction, because the streamwise resolution requirements for streak like motions are less stringent than those for the 
nearfield solution. 

The roughness heights employed in the present set of simulations are small enough that the induced disturbance 
field is purely stationary. Accordingly, the simulations are advanced in time until the last pair of solution snapshots 
separated by a finite time interval can be considered as identical within a specified tolerance. Preliminary 
simulations for the baseline case used a full wavelength domain (-Ao/2 < z < Ao/2), with periodic boundary 
conditions along the spanwise direction. The results of these simulations confirmed the spanwise symmetry of the 
disturbance field, in agreement with the experimental observations of White and Ergin (2003). Accordingly, all of 
the simulations reported in this paper have been carried out by using a computational domain with a spanwise extent 
of -Ao/2 < z < 0 and symmetry conditions at both ends. 

The numerical simulations are initialized with the Blasius similarity solution, shifted appropriately in the wall- 
normal direction in the region of the roughness array. The relevant Blasius profiles are also imposed along the 
inflow boundary, whereas no-slip conditions are used along the entire plate surface. The free-stream boundary 
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conditions enforce an irrotational disturbance behavior near the farfield boundary, and a convective boundary 
condition was imposed along the downstream (i.e., outflow) boundary. The robustness and effectiveness of the 
various boundary conditions were verified in separate computations involving spatial simulations of instability wave 
amplification within a parallel boundary-layer flow (Fig. 2). 

For spectral element calculations, grid convergence can be assessed by reducing the size of the elements and/or 
by increasing the order N of the discretization polynomial within each element. Typically, the latter approach is 
more convenient; this approach was adopted in the present study to evaluate the grid convergence of the computed 
disturbances. In particular, nearfield and farfield solutions obtained with 11th order and 7th order elements, 
respectively, were compared with each other for the case of h = 0.665 mm (i.e., the largest roughness height 
considered in the present study). Figures 3a and 3b illustrate the agreement between streamwise and spanwise 
velocity components sampled from the N=7 and N=ll solutions, respectively, at an x-location near the outflow of 
the nearfield domain. For the most part, a similar agreement was found between the corresponding solutions within 
the farfield domain, except for a noticeable discrepancy near the outflow of the domain. However, the magnitude of 
this discrepancy was small and was confined to the far downstream region wherein the streamwise element size had 
become rather large. Comparison between the nearfield and farfield solutions within the domain of overlap (wherein 
the two domains involve somewhat different grid distributions) provided an additional degree of check with respect 
to the grid convergence of the computed solution, in addition to validating the domain decomposition approach used 
in these simulations. 

Nearly all of the postprocessing of the simulation data was based on interpolating the computed solution to a 
specified output grid involving uniform spacing of points along the spanwise direction. To minimize any loss of 
accuracy, spectral interpolation based on the same degree of polynomial as that used for the simulation itself was 
employed during this interpolation. Fast cosine and sine transforms were used to decompose the interpolated 
solution into Fourier modes; integrated modal energies at the streamwise locations of interest were evaluated using 
the trapezoidal rule for quadrature. 

III. Results for Baseline Configuration 

Despite the relative simplicity of the overall problem, the parametric space involved is rather large, including 
both mean boundary-layer parameters at the roughness location (Re 0 , shape factor H, pressure gradient parameter p, 
etc.) and relevant attributes of the roughness geometry (shape, planform size and aspect ratio, and roughness height, 
i.e. Re k ). Limitations of computing resources precluded a detailed parametric study of roughness-induced transient 
growth. Thus, a limited case study has been performed to validate simulation results against measured data, to 
provide additional data for validating simpler theoretical models such as that of Tumin and Reshotko (2004), and to 
delineate the influence of roughness geometry on the process of transient growth. The idea behind the case study 
was that if simplified theoretical models can be sufficiently improved and validated via the combination of wind 
tunnel experiments and numerical simulations, then those same models can be easily used for a detailed parametric 
study of roughness-induced transient growth. With that in mind, simulation results for the baseline configuration of 
h = 0.57 mm are presented in this section; results for the additional roughness configurations are deferred until 
section IV. Unless stated otherwise, all of the simulation results presented in sections III and IV are based on N = 
11 . 


A. Overall Features 

The salient features of the near-wall flow field in the vicinity of the cylindrical roughness element may be 
gleaned from Fig. 4. The bottom half of the figure displays the surface pressure contours, and the upper half depicts 
the distribution of streamwise velocity component at a small, fixed height just above the surface. By confining the 
range of streamwise velocity contours to negative values, the top half of the figure clearly highlights the two 
separate regions of flow reversal — a primary region behind the roughness element and a smaller one near the 
centerline just in front of the element. The associated pattern of 2D streamlines within this plane is also shown to 
emphasize the direction of the in-plane velocity vectors. 

The surface pressure distribution in the bottom half of Fig. 4 is partially analogous to that predicted by Smith et 
al. (1977) by applying the linearized version of triple-deck theory to an isolated roughness element with a smoother 
(i.e., C 1 ) height distribution. In both their theoretical results and the present simulation, one observes peak pressures 
along the centerline just ahead of the roughness element, lowest pressure above the roughness itself, and a pressure 


5 

American Institute of Aeronautics and Astronautics 



recovery in the wake region. However, the present simulations suggest the pressure recovery to be virtually 
monotonic, i.e., without the small overshoot just behind the roughness element as predicted by the linearized theory. 
The discrepancy between the theoretical predictions of Smith et al. (1977) (as well as Tumin and Reshotko 2004) 
and the simulation results herein may be due to of the following: the differences in roughness element shape and 
spanwise spacing and the effects of nonlinearity which had been neglected in the theory (and which may alter the 
effective roughness shape when the region of flow reversal behind the circular roughness has been taken into 
account). Additional calculations are necessary to help identify the precise cause behind this discrepancy. 

The computed distribution of streamwise shear stress within the wake region confirmed the theoretically 
predicted feature of a streamwise elongated “corridor” involving a pair of high-speed streaks, one on each side of the 
centerline (see Fig. 21 from Tumin and Reshotko 2004). The centerline region displays a low-speed streak 
involving a progressively reducing velocity deficit. Visualization of the 3D flow field suggests that the occurrence 
of high-speed streaks may be attributed to the horseshoe vortex system formed as the spanwise vorticity lines in the 
incoming flow wrap around the cylindrical roughness element in the near-wall region. Velocity contours in the wake 
of the roughness array (Fig. 5) illustrate the longevity of the vortex structures associated with the high-speed streaks. 
The spanwise structure of the wake disturbances is qualitatively similar to the computations of Joslin and Grosch 
(1995) and Woerner et al. (2004) for smoother roughness shapes. 

The streamwise velocity contours at a couple of selected streamwise stations downstream of the roughness array 
are show n in Figs. 6(a) and 6(b). The locations of these stations were chosen to allow a direct comparison with Figs. 
4 and 7, respectively, from White and Ergin (2003). The near-wake contours at x - x 0 = 10 mm (Fig. 6(a)) indicate a 
single dominant low-speed streak near the centerline, similar to that observed in the experiments of Kendall (1982), 
Gaster et al. (1994), and White and Ergin (2003). The high-speed streaks mentioned previously emerge farther 
downstream, as seen from Fig. 6(b). 

Streamwise evolution of the total disturbance energy (i.e., mean square amplitudes integrated across the profile 
at a fixed streamwise location) illustrates the weak decay in spanwise varying perturbations within the far-wake 
region as shown by the solid curve in Fig. 7. To maintain consistency with the data reported by White and Ergin 
(2003), the disturbance energy profiles are integrated with respect to the local similarity variable at each streamwise 
station. This choice of wall-normal length scale masks the effects of boundary-layer growth in the streamwise 
direction. However, the energy of the spanwise varying perturbations is found to decay (albeit more weakly) even 
after the boundary-layer growth has been accounted for in calculating the total disturbance energy at any given 
location. The dashed curve in Fig. 7 is based on a reference state corresponding to the z = -Xq/2 plane, which lies in 
between an adjacent pair of roughness elements. All of the features of the disturbance energy evolution in Fig. 7 are 
in satisfactory agreement with Fig. 1 1 from White and Ergin (2003). 

B. Transient Growth and Comparison with Measured Data 

Having validated the general features of the computed disturbance field against the measured data, we next 
examine the downstream evolution of the individual Fourier harmonics in order to quantify detailed features of the 
transient growth induced by the baseline roughness array with h = 0.57 mm. The disturbance profiles at the 
fundamental spanwise wave number and its first three harmonics at streamwise locations of x - x 0 = 10 mm and x - 
x 0 = 80 mm are shown in Figs. 8(a) and 8(b), respectively. These two locations are the same as those considered in 
Figs. 6(a) and 6(b) discussed previously. A comparison of the modal amplitudes from Figs. 8(a) and 8(b) indicates 
the switchover in the dominant Fourier component from X = Xo at the upstream location to the X = Xq/3 mode farther 
downstream, as previously noted by White and Ergin (2003). The root mean square (rms) disturbance amplitudes at 
both locations agree quite well with the measured data presented in Figs. 6 and 9, respectively, from White and 
Ergin. A scrutiny of the rms profiles at other stations did not reveal any double peaks, supporting the claim by 
White and Ergin that the presence of double peaks in their earlier results (White 2002) was caused by inaccuracies 
related to the wall finding algorithm. 

The modal profiles from White and Ergin (2003) correspond to the power spectral density (PSD) at a given 
spanwise wavelength rather than to the mean-square amplitude over each spectral peak as plotted in Fig. 8. 
Assuming a constant shape and wave number bandwidth across the relevant set of peaks in the measured PSD data, 
the results from Fig. 8 can be compared with the modal profiles in Figs. 6 and 9 of White and Ergin via a simple 
rescaling. Accordingly, we have rescaled the computed mean-square amplitudes in Fig. 8 by a constant factor of 
13.7. This particular scaling factor was chosen to provide an approximate match between the respective peaks in the 
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modal profiles for X = Aoat x - x 0 = 10 mm. Reasonable agreement is obtained between simulation results and the 
measured data, both in regard to the mode shapes and amplitudes for X = Ao, X = Aq/3, and X = Xq/4. A noteworthy 
discrepancy involves what appears to be a greater reduction in the fundamental mode amplitude in the experiment as 
the measurement location is moved from x - x 0 = 10 mm to x - x () = 80 mm. As we observe later, the fundamental 
mode does not decay monotonically between the two locations, which may account for the discrepancy related to the 
decay factor. A comparison between Figs. 8(a) and 8(b) also confirms the outward movement in the peak locations 
corresponding to the modal profiles for X = Aq/3 and X = Xo/4 as one moves downstream from x - x 0 = 10 mm to x - 
x 0 = 80 mm. 


Streamwise evolution of the integrated modal energy for the wave number components X = Xo, X = A</3, and X = 
Xq/4 is shown in Fig. 9(a). All three modes exhibit transient growth before decaying farther downstream; however, 
the fundamental mode undergoes a rapid and substantial decay before the beginning of its algebraic growth (and its 
eventual decay is not captured within the computational domain). The shorter wavelength modes, on the other hand, 
begin to amplify almost immediately in the wake of the array, similar to the experiment. The X = Xq/4 mode reaches 
its peak amplitude before the other modes, and is closely followed by the X = Xq/3 mode which dominates the 
disturbance amplitudes in the intermediate wake region (up to approximately 170 mm downstream of the roughness 
array). Interestingly enough, the minimum of the fundamental mode energy occurs just slightly upstream of the 
maximum in the X = A</3 mode on the scale of this plot. Again, similar to Fig. 7, we have used the local similarity 
length scale to nondimensionalize the wall-normal coordinate r\ in order to permit a visual comparison with the 
measurements of White and Ergin which have been reproduced in Fig. 9(b) for the reader’s convenience. (Note that 
the comparison may be only qualitative at this stage because of the potentially nontrivial conversion from PSD to 
squares of the modal amplitudes.) With the choice of a streamwise varying length scale, the energy metric (E^x)) 
for the fundamental mode appears to reach a plateau near Xe* p = 600 mm (i.e., x - x 0 = 300 mm) even though the 
actual disturbance energy continues to grow throughout the streamwise region included in Fig. 9(a). The continued 
growth of the fundamental mode appears consistent with the predictions of the optimal growth theory ( Andersson et 
al. 1999, Tumin and Reshotko 2001), even though the initial locations considered in those calculations were 
typically much closer to the leading edge than the roughness array considered herein. If the previously published 
optimal growth results were to hold even for the present case, then the fundamental mode X = Xohas the potential to 
amplify until the local spanwise wave number has reached a value of p 0 ~ 0.45, i.e., until approximately x = 1 175 
mm. To determine whether such optimal growth is realized in this case, the present far-field calculations will have 
to be extended upto larger values of x. The additional far-field calculations would also help identify (as well as 
resolve) any connection between a possible lack of grid convergence in the far-downstream portion of the present 
farfield domain and the earlier crossover in the computation between the modal energies for X = Xq and X = Ao/3 as 
the latter mode decays across this region. 

Tumin and Reshotko (2004) made the key observation that a sign reversal of the U k modal profiles in the wake 
of the roughness array creates conditions that locally resemble the theoretically predicted optimal initial profiles 
corresponding to a spanwise array of streamwise vortices. Conventional applications of the optimal growth theory 
consider a streamwise interval that begins at the leading edge (or very close to it), in which case these streamwise 
vortices actually protrude outside of the boundary layer. Tumin and Reshotko presented an optimal growth 
calculation based on a streamwise interval beginning at a finite distance downstream of the leading edge; this 
example case suggests that the optimal initial conditions for transient growth over regions away from the leading 
edge also involve a vortical secondary flow, which is now embedded inside the finite-thickness boundary layer at 
the starting location. The evolution of the modal profiles in Fig. 10 shows that the streamwise velocity perturbations 
at X = Ao indeed undergo a sign change near x - Xo = 35 mm (i.e., x exp = 335 mm), which corresponds to the 
minimum in fundamental mode energy in Fig. 9(a). The secondary flow at this location involves the modal profiles 
shown in Fig. 10(b) and 10(c), which are analogous in shape to those predicted by Tumin and Reshotko. 

Because of the small magnitudes of the wall-normal and spanwise velocity perturbations, it is difficult to 
measure those perturbations in a wind tunnel experiment. However, that difficulty is less severe in the numerical 
simulations; we find that the “initial” disturbance energy based on just the v’ and w’ perturbations at X = Xq\s given 
by J(Vj 2 (Tl)+Wi 2 (n)) dq ~ 5.76x10 6 at x - xo = 36 mm. This translates into an algebraic growth in the modal energy 
Ej(x) by a factor of nearly 191 over the transient growth interval x exp = 336 mm to x exp = 650 mm included in the 
computational domain. Of course, from a physical standpoint, one must keep in mind that the above growth is 
preceded by an even larger disturbance decay upstream of the minimum in Ei(x); hence, the overall significance of 
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transient growth via a spanwise periodic roughness array cannot be determined at this stage. Indeed, a roughness 
array with subcritical heights can influence the transition process in multiple ways. It can modify the stability 
characteristics via both strong but localized disturbances in the vicinity of the roughness array and weaker but 
persistent disturbances in the transient growth region (see Cossu and Brandt (2002) for a theoretical study of the 
latter mechanism). Furthermore, it is well-known that the roughness array will also serve as a catalyst for the 
generation of both 2D and 3D Tollmien-Schlichting instability waves, which may or may not begin to grow until 
much later depending on the combination of roughness location and spacing and the frequencies of unsteady free- 
stream disturbances. 

IV, Influence of Roughness Height, Diameter, and Shape on Transient Growth 

Additional simulations were performed to discern the effects of roughness geometry on the transient growth 
phenomenon. A brief outline of the pertinent findings from these simulations is provided in this section. 

A. Effect of Roughness Height 

According to White and Ergin’s data (2003), both the transient growth characteristics and the subsequent 
exponential decay of the X - A</3 mode are approximately uniform across the range of roughness heights considered 
in their measurements. Yet, the modal amplitude at any fixed station downstream of the roughness array is a 
nonlinear function of the roughness height parameter. Both of these findings are qualitatively consistent with 
previous theoretical development for the closely related problem of roughness-induced excitation of inviscid, 
stationary vortex instabilities in a boundary layer (Choudhari and Duck, 1996). The asymptotic theory showed that 
the nearfield dynamics for the roughness elements under consideration becomes nonlinear at relatively small 
roughness heights compared with the boundary layer thickness, but that the disturbance dynamics farther 
downstream of the roughness elements remains linear to the leading order (until, of course, the exponential growth 
of those inviscid instabilities causes nonlinear effects to come into play once again). In other words, the influence of 
nonlinearity is confined to the output of the receptivity process, i.e., to determining the scaling of the effective initial 
amplitude as a function of the roughness height, and that the initial growth of the generated vortex/streak mode is 
independent of the roughness height parameter to the leading order. 

White and Ergin were able to correlate their data for the transient growth of X = XqB and X = Xq/4 modes with a 
quadratic fit of the form E 3 , E 4 <*= Re k 2 at any fixed x. To assess this scaling behavior within the computational 
framework, additional simulations were carried out for h = 0.665 mm (i.e., Re k * 162). A comparison of the 
nearfield solutions for Re k ~ 162 and the baseline case (Re k -119) presented in section III showed that the 16.7% 
increase in roughness height resulted in approximately 24% percent increase in the length of the reversed flow 
region along z=0 behind the roughness element. Modal analysis of the wake disturbances (Fig. 11) suggests that the 
Re k 2 scaling applies to modal energies of all four modes under consideration (i.e., X = Ao, Xq/2, A</3, and Xq/4). 

White and Ergin (2003) did not provide any rationale for choosing the particular scaling noted above. The Re k 2 
scaling for the modal energy translates into an amplitude scaling of the form A TGM ~ C Re k , where A TGM is a 
measure of the streamwise velocity perturbation associated with a given mode undergoing transient growth (e.g., 
Ajcm could be interpreted as max(|f/*(T])|) within the region of transient growth). This amplitude scaling is 
equivalent to A TG m ~ C 2 h 2 in view of the approximately linear behavior of the Blasius profile over the range of 
roughness heights included in the correlation. Of course, the purely quadratic scaling for the perturbation 
amplitudes cannot hold at sufficiently small roughness heights; therefore, the truncated form of regular perturbation 
expansion, A TGM ~ Cj h + C 2 h 2 , must be considered in general. Numerical simulations are well-suited to 
characterize the linear behavior at small heights (wherein the experimental measurements might encounter 
difficulties in maintaining an adequate signal-to-noise ratio because of the reduced perturbation amplitudes) as well 
as to help delineate the transition from a linear behavior to the quadratic dependence observed in experiments and 
simulations when the Re k parameter is sufficiently large. The parameter study required for this characterization is 
currently under way. 

B. Effect of Ratio of Roughness Diameter to Array Spacing 

To help understand the factors responsible for the absence of transient growth in the X = Xq/2 mode, we carried 
out an additional simulation for the larger roughness radius of r 0 = Xq/4 (i.e., roughness diameter of one-half the 
array spacing). To avoid the complexity of the split domain calculation, the simulation for the larger planform-size 
roughness array was based on a single structured grid of elements within a rectangular box volume. Also, the no- 
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slip conditions over the roughness array were transferred to the plate surface y = 0 by using a perturbation expansion 
based on a linear Taylor series expansion in the roughness height parameter. Since the solver is based on the full 
Navier-Stokes equations, the height of the roughness element was reduced to h = 0.12 mm in order to reduce the 
influence of the nonlinear terms yet avoid a loss of accuracy in analyzing the perturbations. A half-wavelength 
cosine fairing was used to round off the edges of the cylindrical roughness element over the fringe region r/r 0 = 1± 
0.25. Based on the grid convergence study for the baseline configuration, the discretization order was also reduced 
from N = 1 1 to N = 8 for this configuration; however, no additional grid convergence checks were performed. 

Because the spacing between the high-speed streaks scales with the roughness diameter, the X = Xo/2 mode 
emerges in this case as the dominant Fourier harmonic in a large portion of the wake region. However, similar to 
the baseline configuration, this mode decays monotonically, without any evidence of transient growth. The X = Xq 
and X - Ao/3 modes do show an algebraic amplification; however, their amplitudes are significantly weaker than the 
baseline configuration (r 0 = Xq/6 ). Furthermore, the growth of the fundamental mode begins relatively farther 
downstream, and (unlike the baseline configuration) the X = Xo/4 mode does not grow at all. These findings indicate 
the significance of roughness diameter (in relation to array spacing) as a significant parameter in determining the 
quantitative aspects of roughness-induced transient growth. 

C. Roughness Elements with Square Planform and Comparison with Theory of Tumin and Reshotko (2004) 

* 

The theoretical predictions of Tumin and Reshotko (2004) based on a linear model are in qualitative agreement 
with the measurements of White and Ergin (2003) and the simulation results despite the obviously nonlinear regime 
of roughness heights in both experiment and simulations. However, there is a notable discrepancy between theory 
and measured data with regard to the onset of transient growth in the X = Xq/ 3 and X = Xq/4 modes. The algebraic 
growth of these modes begins almost immediately behind the roughness array in the experiment, whereas it is 
postponed until approximately 150 mm downstream of the roughness array according to the theoretical predictions. 
To help ascertain whether the discrepancy is related to differences in planform shapes of the roughness elements (a 
square shape in the calculations of Tumin and Reshotko, versus circular disks in the experiment), an additional 
simulation was carried out for square shaped elements with the same width as the roughness diameter in the baseline 
configuration (i.e., r 0 = Ao/6). In an attempt to simulate the linear behavior modeled by theory, the modified 
simulation procedure based on a transfer of roughness boundary condition to the plate surface, as described in 
section IV-B, was also applied to this particular configuration. The results indicated that the transient growth 
features for square elements were analogous to those for the circular elements of comparable size, in that the X = 
Ao/3 and X = X</4 modes again began to amplify just downstream of the roughness array. Additional work is 
therefore necessary to identify the cause behind this discrepancy. However, once the theoretical model has been 
validated (and refined if necessary) for a simple configuration of this type, such models could be extended to the 
difficult case of distributed, sand-grain roughness (which will require rather prohibitive resources for numerical 
simulations) as well as being used as a predictive tool for roughness-induced transient growth. 

V. Concluding Remarks 

Direct numerical simulations have been carried out to shed further light on the phenomenon of roughness- 
induced transient growth induced by a spanwise periodic array of roughness elements. For roughness arrays with 
Re k =119 and Re k « 162, the computed features of transient growth are consistent with those measured in the 
experiment of White and Ergin (2003), including the nonlinear scaling of perturbation energy with the roughness 
height parameter. The agreement between simulation results and the measured data also confirms the accuracy of the 
spanwise-global wall-finding algorithm used by White and Ergin (2003) over the local wall-finding procedure used 
previously by White (2002). Good quantitative agreement has been shown between the computed and the measured 
rms disturbance amplitudes for the baseline configuration of Re k = 119. The agreement pertaining to modal 
amplitudes is somewhat less satisfactory, and additional work is required to help identify the factors underlying this 
discrepancy. Similarities between the transient growth due to roughness array composed of square shaped and 
circular roughness elements indicate the likely significance of roughness size relative to array spacing as a primary 
control variable, at least for the class of cylindrical roughness elements examined herein. 

A number of relevant issues still need to be addressed, such as in-depth comparison with both experiments and 
linearized theory, calculations for additional roughness heights so as to better quantify the overall scaling of 
disturbance amplitudes (and, secondarily, the magnitudes of transient growth) with Re k . Further calculations would 
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also help delineate the critical range of roughness heights beyond which the flow becomes unsteady and vortex 
shedding occurs as the precursor to the formation of turbulent spots (Acarlar and Smith 1987, Klebanoff et al. 1992). 
An additional, relatively straightforward extension of the present simulations would involve the roles of planform 
aspect ratio (i.e. elliptical instead of circular disk roughness elements) and spanwise asymmetry of the planform 
shape of roughness elements. 

While the work provides an independent confirmation for the occurrence of transient growth due to surface 
roughness, it also underscores the weakness of the peak disturbance amplitudes involved in the transient growth 
process. Although it seems unlikely that the role of nearfield disturbances can be diminished via proper selection of 
roughness geometry, investigations focused on the identification of optimal roughness distributions that maximize 
the amplitude of the wake disturbances while keeping the amplitudes of the nearfield disturbances to a minimum 
would also be beneficial. 

Finally, a number of recent studies (Bakchinov et al. 1995, Asai et al. 2002, Cossu and Brandt 2002, Wu and 
Luo 2003) have addressed the influence of stationary streaks and optimal growth disturbances on the stability of the 
unperturbed boundary layer and/or the interaction between stationary streaks and the unsteady Tollmien-Schlichting 
instabilities in the boundary layer. Additional studies of this type in the specific context of roughness-induced 
transient growth would be helpful, especially from the standpoint of placing the roughness-induced transient growth 
in a correct physical perspective. 
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(a) Schematic of the combined near- and far-zone grids. 


(b) Close-up of the grid near the roughness 
element, including a symmetry-plane cut. 


Figure 1. Spectral element mesh. 






Figure 2. Validation of spectral element code for Tollmien-Schlichting wave propagation (F=70xl0' 6 ). 

Top: contours of v’ downstream of source (an unsteady suction-blowing strip at Re x = 0.5x1 0 6 ). 
Bottom: Magnitude of v’ vs. x at y=l, compared with linear stability theory (dashed line). 
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(a) u profiles (b) w profiles 

Figure 3. Grid convergence of nearfield solution for largest height roughness array (h = 0.665 mm): velocity profiles 
at x - xq = 0.008 m. Velocity components u and w (measured in m/s) are plotted as functions of span wise 
coordinate z (in meters) for various heights above the surface. Solid lines denote N = 1 1 solution, whereas 
dashed lines indicate the N = 7 solution. 
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Figure 4. Flow features near roughness array: Figure 5. Stationary vortex structures in the wake of spanwise 

surface pressure distribution (lower half) and periodic roughness array (baseline configuration). 

streamwise shear stress (upper half). Warmer 

shades denote peak pressures or large positive 

shear stress; cooler shades denote the opposite 

conditions. 



(a) x-x 0 = 10 mm (a) x-x 0 = 80 mm 

Figure 6. Streamwise velocity contours at x-x 0 = 10 mm and 80 mm. respectively (contour labels indicate velocity in m/s). 
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Figure 7. Evolution of total disturbance energy behind baseline roughness array (h = 0.57 mm). 




Figure 8. RMS and modal disturbance profiles at x-xO = 10 mm and 80 mm, respectively. The modal profiles correspond 
to square of the peak amplitude in the corresponding Fourier component, scaled by a constant factor of 13.7. 



(a) Modal energies ju k 2 (q)dr|) based on simulation (a) Measured data from White and Ergin (2003) 

Figure 9. Streamwise evolution of modal energy in dominant Fourier components. 
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Figure 1 0. Modal profiles for the fundamental spanwise mode at selected streamwise locations behind roughness array 
(h = 0.57 mm). 



Figure 1 1 . Effect of roughness height on transient growth: normalized modal energies E k (x)/Re k 2 (k=l-4), plotted 
as functions of x for baseline array (Re k =119) and the larger-height array with Re k = 162. 
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